R for Data Science

library(tidyverse)
library(plotly)
library(ggpubr)
library(rstatix)

Attaching package: ‘rstatix’

The following object is masked from ‘package:stats’:

    filter

Read in the gapminder_clean.csv data as a tibble using read_csv.

gapminder_clean <- read_csv("gapminder_clean.csv")

New names:
* `` -> ...1
Rows: 2607 Columns: 20
── Column specification ────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Delimiter: ","
chr  (2): Country Name, continent
dbl (18): ...1, Year, Agriculture, value added (% of GDP), CO2 emissions (metric tons per capita), Domestic credit provided by financial sector ...

ℹ Use `spec()` to retrieve the full column specification for this data.
ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

Filter the data to include only rows where Year is 1962 and then make a scatter plot comparing ‘CO2 emissions (metric tons per capita)’ and gdpPercap for the filtered data.

gapminder_1962 <- filter(gapminder_clean, Year==1962)

ggplot(gapminder_1962, aes(x=`CO2 emissions (metric tons per capita)`, y=gdpPercap)) + geom_point()

On the filtered data, calculate the correlation of ‘CO2 emissions (metric tons per capita)’ and gdpPercap. What is the correlation and associated p value?

```r
cor.test(gapminder_1962$`CO2 emissions (metric tons per capita)`, gapminder_1962$gdpPercap)

<!-- rnb-source-end -->

<!-- rnb-output-begin eyJkYXRhIjoiXG5cdFBlYXJzb24ncyBwcm9kdWN0LW1vbWVudCBjb3JyZWxhdGlvblxuXG5kYXRhOiAgZ2FwbWluZGVyXzE5NjIkYENPMiBlbWlzc2lvbnMgKG1ldHJpYyB0b25zIHBlciBjYXBpdGEpYCBhbmQgZ2FwbWluZGVyXzE5NjIkZ2RwUGVyY2FwXG50ID0gMjUuMjY5LCBkZiA9IDEwNiwgcC12YWx1ZSA8IDIuMmUtMTZcbmFsdGVybmF0aXZlIGh5cG90aGVzaXM6IHRydWUgY29ycmVsYXRpb24gaXMgbm90IGVxdWFsIHRvIDBcbjk1IHBlcmNlbnQgY29uZmlkZW5jZSBpbnRlcnZhbDpcbiAwLjg5MzQ2OTcgMC45NDg5NzkyXG5zYW1wbGUgZXN0aW1hdGVzOlxuICAgICAgY29yIFxuMC45MjYwODE3IFxuIn0= -->

Pearson's product-moment correlation

data: gapminder_1962\(`CO2 emissions (metric tons per capita)` and gapminder_1962\)gdpPercap t = 25.269, df = 106, p-value < 2.2e-16 alternative hypothesis: true correlation is not equal to 0 95 percent confidence interval: 0.8934697 0.9489792 sample estimates: cor 0.9260817

<!-- rnb-output-end -->

<!-- rnb-chunk-end -->


<!-- rnb-text-begin -->

There is a positive  correlation of 0.926 (p-val <2.2e-16).

**On the unfiltered data, answer "In what year is the correlation between 'CO2 emissions (metric tons per capita)' and gdpPercap the strongest?" Filter the dataset to that year for the next step...**


<!-- rnb-text-end -->


<!-- rnb-chunk-begin -->


<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxubWF4X3llYXJfY29yIDwtIGdhcG1pbmRlcl9jbGVhbiAlPiUgZ3JvdXBfYnkoWWVhciklPiUgXG4gIHN1bW1hcml6ZSh5ZWFyX2Nvcj0gY29yKGBDTzIgZW1pc3Npb25zIChtZXRyaWMgdG9ucyBwZXIgY2FwaXRhKWAsIGdkcFBlcmNhcCwgIHVzZT1cImNvbXBsZXRlLm9ic1wiKSkgICU+JSB0b3BfbihuPTEsIHd0PXllYXJfY29yKSAlPiUgcHJpbnQoKVxuYGBgIn0= -->

```r
max_year_cor <- gapminder_clean %>% group_by(Year)%>% 
  summarize(year_cor= cor(`CO2 emissions (metric tons per capita)`, gdpPercap,  use="complete.obs"))  %>% top_n(n=1, wt=year_cor) %>% print()

gapminder_1967 <- filter(gapminder_clean, Year==1967)

The correlation is strongest in 1967 (0.9388).

Using plotly, create an interactive scatter plot comparing ‘CO2 emissions (metric tons per capita)’ and gdpPercap, where the point size is determined by pop (population) and the color is determined by the continent. You can easily convert any ggplot plot to a plotly plot using the ggplotly() command.

p <- ggplot(gapminder_1967, aes(x=`CO2 emissions (metric tons per capita)`, y=gdpPercap, size=pop, color=continent)) + geom_point() 

ggplotly(p)

What is the relationship between continent and ‘Energy use (kg of oil equivalent per capita)’? (stats test needed)

gapminder_energy <- select(gapminder_clean, continent, `Energy use (kg of oil equivalent per capita)`) %>% filter(complete.cases(.))
norm_test_energy <- gapminder_energy %>% group_by(continent) %>% summarize(shap_test_p_val= shapiro.test(`Energy use (kg of oil equivalent per capita)`)$p.value)
res.kruskal <-kruskal_test(`Energy use (kg of oil equivalent per capita)` ~ continent, data = gapminder_energy)
pwc <- gapminder_energy %>% 
  wilcox_test(`Energy use (kg of oil equivalent per capita)` ~ continent, p.adjust.method = "bonferroni")

pwc <- pwc %>% add_xy_position(x = "continent")
ggboxplot(gapminder_energy, x = "continent", y = "Energy use (kg of oil equivalent per capita)") +
  stat_pvalue_manual(pwc, hide.ns = TRUE) +
  labs(
    subtitle = get_test_label(res.kruskal, detailed = TRUE),
    caption = get_pwc_label(pwc)
    )

# ggplot(gapminder_energy, aes(x=continent, y=`Energy use (kg of oil equivalent per capita)`)) + geom_boxplot() + stat_compare_means(method="kruskal.test")

Data is not normally distributed

Is there a significant difference between Europe and Asia with respect to ‘Imports of goods and services (% of GDP)’ in the years after 1990? (stats test needed)

gapminder_af_1990 <- filter(gapminder_clean, Year >1990 & (continent=="Europe" | continent =="Asia"))
norm_test_1990 <- gapminder_af_1990 %>% group_by(continent) %>% summarize(shap_test_p_val= shapiro.test(`Imports of goods and services (% of GDP)`)$p.value)
norm_test_1990 

The normality test shows that the data are not normally distributed. Therefore, to compare both groups we need a non-parametric test (like the Mann–Whitney or Wilcoxon rank-sum test).

ggplot(gapminder_af_1990, aes(x=continent, y=`Imports of goods and services (% of GDP)`, fill=continent)) + geom_violin() + stat_compare_means(method = "wilcox.test", label.x = 1.4, label.y = 200)+
  theme(legend.position="none")

No, there is not a significant difference between Europe and Asia with respect to ‘Imports of goods and services (% of GDP)’ in the years after 1990.

What is the country (or countries) that has the highest ‘Population density (people per sq. km of land area)’ across all years? (i.e., which country has the highest average ranking in this category across each time point in the dataset?)

```r
count_high_pop_dens <- gapminder_clean %>% group_by(`Country Name`) %>% summarize(av_pop_dens = mean(`Population density (people per sq. km of land area)`, na.rm = T)) %>% top_n(n=1, wt=av_pop_dens) %>% print()

<!-- rnb-source-end -->

<!-- rnb-frame-begin eyJtZXRhZGF0YSI6eyJjbGFzc2VzIjpbInRibF9kZiIsInRibCIsImRhdGEuZnJhbWUiXSwibnJvdyI6MSwibmNvbCI6Miwic3VtbWFyeSI6eyJBIHRpYmJsZSI6IjEgw5cgMiJ9fSwicmRmIjoiVWtSWU13cFlDZ0FBQUFNQUF3WURBQU1GQUFBQUFBVlZWRVl0T0FBQUJBSUFBQUFCQUFRQUNRQUFBQUY0QUFBREV3QUFBQUlBQUFBUUFBQUFBUUFFQUFrQUFBQVFUV0ZqWVc4Z1UwRlNMQ0JEYUdsdVlRQUFBQTRBQUFBQlFNekdCSUhIZ0YwQUFBUUNBQUFBQVFBRUFBa0FBQUFGYm1GdFpYTUFBQUFRQUFBQUFnQUVBQWtBQUFBTVEyOTFiblJ5ZVNCT1lXMWxBQVFBQ1FBQUFBdGhkbDl3YjNCZlpHVnVjd0FBQkFJQUFBQUJBQVFBQ1FBQUFBbHliM2N1Ym1GdFpYTUFBQUFOQUFBQUFvQUFBQUQvLy8vL0FBQUVBZ0FBQUFFQUJBQUpBQUFBQldOc1lYTnpBQUFBRUFBQUFBRUFCQUFKQUFBQUNtUmhkR0V1Wm5KaGJXVUFBQUQrQUFBRUFnQUFBQUVBQkFBSkFBQUFCMjl3ZEdsdmJuTUFBQUlUQUFBQUJRQUFBQkFBQUFBQkFBUUFDUUFBQUFGeUFBQUFFQUFBQUFFQUJBQUpBQUFBRDNWdWJtRnRaV1F0WTJoMWJtc3RNUUFBQUFvQUFBQUJBQUFBQUFBQUFBMEFBQUFCQUFBQUFRQUFBQTBBQUFBQkFBQUFBZ0FBQkFJQUFBTC9BQUFBRUFBQUFBVUFCQUFKQUFBQUJtVnVaMmx1WlFBRUFBa0FBQUFGYkdGaVpXd0FCQUFKQUFBQURuSnZkMjVoYldWekxuQnlhVzUwQUFRQUNRQUFBQXB5YjNkekxuUnZkR0ZzQUFRQUNRQUFBQXBqYjJ4ekxuUnZkR0ZzQUFBQS9nQUFBUDQ9In0= -->

<div data-pagedtable="false">
  <script data-pagedtable-source type="application/json">
{"columns":[{"label":["Country Name"],"name":[1],"type":["chr"],"align":["left"]},{"label":["av_pop_dens"],"name":[2],"type":["dbl"],"align":["right"]}],"data":[{"1":"Macao SAR, China","2":"14732.04"}],"options":{"columns":{"min":{},"max":[10],"total":[2]},"rows":{"min":[10],"max":[10],"total":[1]},"pages":{}}}
  </script>
</div>

<!-- rnb-frame-end -->

<!-- rnb-chunk-end -->


<!-- rnb-text-begin -->


The country with the highest 'Population density (people per sq. km of land area)' across all years in Macao SAR, China.

**What country (or countries) has shown the greatest increase in 'Life expectancy at birth, total (years)' since 1962?**


<!-- rnb-text-end -->


<!-- rnb-chunk-begin -->


<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuY291bnRfaW5jcl9saWZlX2V4cGVjdCA8LSBnYXBtaW5kZXJfY2xlYW4gJT4lIHNlbGVjdChgQ291bnRyeSBOYW1lYCwgYExpZmUgZXhwZWN0YW5jeSBhdCBiaXJ0aCwgdG90YWwgKHllYXJzKWAsIFwiWWVhclwiKSAlPiUgZmlsdGVyIChZZWFyPT0xOTYyIHwgWWVhcj09IG1heChZZWFyKSkgJT4lIGdyb3VwX2J5KGBDb3VudHJ5IE5hbWVgKSU+JSBzdW1tYXJpemUoaW5jX2xpZmVfZXhwZWN0ID0gYExpZmUgZXhwZWN0YW5jeSBhdCBiaXJ0aCwgdG90YWwgKHllYXJzKWBbWWVhcj09bWF4KFllYXIpXS1gTGlmZSBleHBlY3RhbmN5IGF0IGJpcnRoLCB0b3RhbCAoeWVhcnMpYFtZZWFyPT1taW4oWWVhcildKSAlPiUgdG9wX24obj0xLCB3dD1pbmNfbGlmZV9leHBlY3QpICU+JSBwcmludCgpXG5gYGAifQ== -->

```r
count_incr_life_expect <- gapminder_clean %>% select(`Country Name`, `Life expectancy at birth, total (years)`, "Year") %>% filter (Year==1962 | Year== max(Year)) %>% group_by(`Country Name`)%>% summarize(inc_life_expect = `Life expectancy at birth, total (years)`[Year==max(Year)]-`Life expectancy at birth, total (years)`[Year==min(Year)]) %>% top_n(n=1, wt=inc_life_expect) %>% print()

Maldives is the country that shows a greatest increase in life expectancy at birth since 1962.